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Abstract 

The properties of the ground state of liquid "^He are studied using a 
correlated basis function of the form Y[i<j '^(''^ij) ■ Here, iplr) is chosen 
as the exact solution of the Schrodinger equation for two ^He atoms. A 
hard-sphere plus an attractive square well is used as the interaction po- 
tential between ^He atoms. The pair distribution function is calculated 
using approximate integral methods, namely the Percus-Yevick (PY) 
equation and Hypernetted Chain (HNC) approximation. The values 
thus obtained are used to calculate the ground state energy, which is 
found to be —4.886 K using the PY equation. The liquid structure fac- 
tor is also obtained using the pair distribution function. The values for 
the pair distribution function and liquid structure factor are compared 
with experimental results and earlier theoretical calculations. 

1. Introduction 

The ground state properties of rare gas fluids have been long studied (P - [3]). 
Many of these properties are determined by the knowledge of pair distribution 
function, g{r), which can be measured experimentally. There have been 
many efforts to calculate the pair distribution function over the years. Most 
of the theoretical efforts employed methods mostly used in classical physics, 
such as Monte Carlo ([1] - [7|). Alternately, the pair distribution function is 
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determined by solving a nonlinear integral equation such as Percus-Yevick 
(PY) |H] or the Hypernetted Chain approximation (HNC) [H] . Most of these 
calculations are done in the classical physics domain. 

For liquid helium, quantum effects are important and has to be incorpo- 
rated in the theory. Some of the calculations for obtaining the pair distribu- 
tion function of liquid helium employed Quantum Monte Carlo methods ( [TU] 
- [Hj) but most of the pair distribution function studies used only classical 
methods ([15] - [IZ]). The properties of liquid helium can be understood to a 
great extent in terms of short-range correlations emanating from the strong 
interactions between particles. The interactions are usually represented by 
a central pair potential of short range, the most common form being the 
Lennard-Jones(LJ) potential. An alternate to the LJ potential is the hard- 
sphere potential. 

A major success of the hard-sphere potential was that the exact solution 
of the PY equation can be obtained for this potential ([18] - [20]). Several 
later works have extended the calculation by either modifying the hard-sphere 
potential (e.g. considering an extended hard core potential [2T]) or adding 
another form of potential outside the hard-sphere (e.g. an attractive square- 
well potential). These calculations have been classical ([22] - [25]). 

Here we present a quantum mechanical solution to g{r) for a fluid with 
the hard-sphere potential and an attractive square well outside of it. The 
short-range correlations in liquid helium are commonly treated by writing 
the many-body wavefunction as a product of pair functions, known as the 
Jastrow function [26]. We solve exactly the two-body Schrdinger equation 
for the hard-sphere plus attractive square well potential, and then use the 
pair wave function to construct the many-body Jastrow wavefunction. The 
important advantage of such a wavefunction is the formal analogy between its 
energy expectation values with the configuration-space integrals encountered 
in classical equilibrium statistical mechanics. We calculate pair distribution 
function, g{r) of liquid ^He using two approximate integral equation meth- 
ods, the Percus-Yevick (PY) and Hypernetted Chain (HNC) approximation. 
We use the values of the pair distribution function to obtain the ground 
state energy and liquid structure factor, S{k). We compare our theoretically 
obtained values of g{r) and S{k) to the experimental results ([27] - [29] ) 
and also to earlier classical calculations ([15] - [IZ]). Our theoretical results 
qualitatively matches the experimental data. However, one needs to use a 
more accurate form for the potential, as for example, the Aziz potential [30] 
to obtain a better quantitative agreement. One logical and very interesting 
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extension of this work would be to investigate if the wavefunction obtained 
can describe supersohdity. 

2. Theory 

2.1 Ground-State Wavefunction 

The Hamiltonian for a system of ^He atoms of mass m in a volume Q, 
interacting through a central pair potential V^(|r|) is given by, 

1=1 i<j 

The best wavefunction for liquid ^He is the correlated basis function of the 
Jastrow form. For a many-particle fluid, the wavefunction can be written as 

N 
i<j 

= exp^ln[V^(|fi -fjl)] (3) 

i<j 

Now, we must choose a reasonable form for the pair function ■0(r); this func- 
tion should be small for short distances and should approach a constant for 
large distances. At small distances, where the two particles interact strongly 
the pair function is not expected to be very different from the solution of the 
two-body problem. If the potential function between pairs of ^He atoms has 
the form of a hard-sphere with an attractive square well, 

{oo r < a 
-Vo a<r<b (4) 
b<r , 

then the Schrodinger equation for the two atom wavefunction, in relative 
coordinates, is 

EtP{r) = [ — + V{r)]^P{r). (5) 



3 



The eigenfunctions (/ = 0) can be written as 

{0 r < a 

Asm[p{r — a)] a < r < b (6) 
sin(A;r + 5) b < r, 

p = k + — -, E = . (7) 

Now, two ^He atoms form a very weakly bound state, with E ^ 0. Hence, 
we can approximate the eigenfunction as 

{0 r < a 

siB.[po{r — a)] a < r < b (8) 
1 b <r, 

Pi - ^. (9) 

TT 

sm\po{b-a)] = l,po{b-a) = -. (10) 

2.2. Integral Equation Methods 

The pair distribution function is defined as 

^(^1, ^2) = [YI SiRi-ri)S{R2-fj)^^dn ...dfN / ^^dn . . . dTN. (11) 

Since, we have a tranlslationally invariant system, the pair distribution func- 
tion is a function of the relative coordinate, (i?i — R2) only. In addition, 
g{Ri — R2) is spherically symmetric in the liquid state, hence, g{Ri — R2) = 
g{\Ri — -R2I) = g{f)- The pair distribution function is normalized in such a 
way that g(r) = 1 for large r and it satisfies the sum rule, 

n j [l-g{r)\dr = 1. (12) 

One important advantage of using correlated basis functions is that the form 
of the diagonal density matrix is mathematically identical to that of 
a classical fluid. The classical problem has been extensively studied, and 
accurate computational methods are available which work well for liquids 
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of neutral atoms. The most successful methods are based on the Percus- 
Yevick(PY) integral equation and the Hypernetted Chain(HNC) approxima- 
tion. The FY equation for the pair distribution function of a system of ^He 
atoms with pair wavefunction described by EqlH] is 



1+nj dy[g{y) - l]g{r- y)[l - ^,^_!_ _,^ ]. (13) 



Performing the integration over angular variables and defining = g'{r), 
EqJT3] can be written as 

2TTn /"~ r~^y 
g'(r) = 1 + / ydy[ij\y)g'iy) - 1] / zdzg'iz)[^\z) - 1]. (14) 

^ Jo J\r-y\ 

On the other hand, the Hypernetted Chain integral equation (HNC) for the 
pair distribution function of the system of interest is given by 

log^ = n| dy[g{y) " ^Mr- v) - I - hgj^^]. (15) 

Using a similar definition for g'{r) as in FY equation and after performing 
the integration over angular variables, we have 

logg'ir) = / ydy[ij\y)g'iy) - 1] / dzz[tP'iz)g'{z) - 1 - log^'(^)] 

'0 J\r-y\ 

(16) 
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2.3. Ground State Properties 

The ground-state energy Eq is given by the expectation values of the Hamil- 
tonian in EqUl 



Eo = {H) = J ^H^dn . . . dr^ j j "^^dri ...dr^. 



(17) 



If the many-body wavefunction, has the form shown in Eqj2], it is straight- 
forward to show that 

I ^H^idfi ...drN= — In ipivij) + V{rij)] "^^dn ...df^ (18) 



In terms of the pair distribution function, the potential energy per particle 
can be written as 

{P.E.) _ n 

N ~ 

and in our case, it is reduced to. 



- j dhV{r)g{r), (19) 



{P.E.) 



N 

The kinetic energy per particle can be obtained using 



27rnK) / drr'^g{r). (20) 

J a 



{K.E.) _ nh' /• 3 d{\n[i;{r)]) dgjr) 



N 4m ./ dr dr 

and, in our case, EqI2I] becomes 

- — ——^ = — — (po [ drg'{r)r^ — 2po j dr g\r)r si\i[pQ{r — a)] cos[pQ{r — a)]) . 

* ^ J a J a 

(22) 

The liquid structure function S{k) (for 7^ 0) is related to the pair distribu- 
tion function g{r) by the following relation 

S{k) = 1 + n j dr[g{r) -l]exp{-ik.r). (23) 

After carrying out the angular integration, the expression for liquid structure 
factor becomes 
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Attti 

S{k) = 1 + ^/ rdr[g{r)-l]sm{\k\\f\). (24) 
|A;| 



3. Method 

The PY integral equation (Eq JT^ and the HNC approximation integral equa- 
tion (Eq IT6|) are solved self-consistently for r < i? = 100 A. The equilib- 
rium density of the system of ^He atoms is chosen to be n = 0.0218 A or 
2.18 X 10^^ atoms/cm^. The hard-sphere radius for the potential function 
is taken to be a = 2.6 A, the width and depth are chosen as 6 = 4 A and 
Vq = —15.26 K respectively. 
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In any numerical method it is necessary to truncate the infinite integrals 
at some stage. These integrals are replaced by a finite sum, assuming that 
g{r) ^ 1 beyond some large but finite radial distance R. We carried out two 
calculations for the self-consistent PY equation (Eq lT^ . one using i? = 30 A 
and another i? = 100 A to investigate the effect of the size of the cutoff. The 
maximum difference between the two sets of g{r) values for all distances is 
3e~^. For the HNC equations, we performed three calculations using R = 50 
A, 60 A and 70 A. The three sets of g{r) values differ not more than by 0.005 
for all distances. We used Simpson's three-point approximation to evaluate 
the integrals at every step of the iteration, using N points evenly distributed 
over the range {0,R). The step size chosen for the Simpson grid was 0.01, 
giving the value of as 10000 for R = 100 A. To investigate the influence of 
the step size on the accuracy of the results we performed two calculations for 
i? = 20 A, one using 2000 points and one using 4000 points. It is found that 
the two sets differ nowhere by more than 0.002. The starting value of g'{r) 
is chosen to be 1 for all values of r. The input values for the next iteration 
is calculated according to the mixing formula 

g'^ir) = a9t'\r) + (1 - c^gf^ir^. (25) 

where a = 0.1 is used to achieve desired convergence. The iterations are 
assumed to give convergence when value of the residual Res is less than 1, 
where Res deflned by 

Res = Y.[g'!2,-gt;'^]\ (26) 

4. Results 

The values of the pair distribution function for ^He obtained using PY 
(Eq jMl) and HNC equations (Eq IT6!) are shown in Fig{Tl The pair distri- 
bution function exhibits familiar features, e.g. g{r) — > 1 for large distances. 
The results satisfy the sum rule within numerical accuracy. Values obtained 
using the PY equation yield a sum of 0.939701, while the values obtained 
using HNC equation give 0.811651. Figure [T] shows our results in comparison 
with earlier classical calcultaions ([I5]-[T7]) as also the x-ray [27] and neutron- 
diffraction data [28l [29]. The different symbols represent the experimental 
data while the different lines represent theoretical calculations. 
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We would like to point out here that while the earlier classical calculations 
yield values which are lower than the experimental data, our values for the 
pair distribution function lies above the experimental results. We would 
ascribe this discrepancy to the choice of parameters for the potential function. 
We used fixed values for the parameters a, b, Vb throughout our calculation. 
Some of the earlier papers (classical calculations) [121 EI] have discussed the 
effect of varying the parameters to both the pair distribution function and 
the liquid structure factor. The earlier paper [16] has shown that the peak of 
the pair distribution function increases with the increase in the hard-sphere 
radius. And the recent article [31] has shown that the peak value of the 
liquid structure factor increases with the increase of width and depth of the 
square well part of the potential. Hence, we need to choose the values of 
the parameter carefully in order to obtain a good quantitative agreement 
between the theoretically obtained pair distribution function values and the 
experimental data. The better way would be to obtain the parameters using 
variational calculation. 

The numerically obtained values of the pair distribution function is uti- 
lized to obtain the ground state energy for ^He. The average value for the 
potential energy per particle and the average value of the kinetic energy per 
particle are —41.67 K and 36.78 K respectively. Hence, the average ground 
state energy per particle is found to be —4.89 K for the pair distribution func- 
tion obtained using PY equation. This value differ from the experimentally 
obtained ground state energy of ^He atoms by about 31%. 

Using our results of the pair distribution function and the relation be- 
tween the pair distribution function and the hquid structure factor (Eq j23|l . 
we now calculate the liquid structure factor of ^He. The results are shown 
in FigJ2] along with earlier theoretical calculations as well as the x-ray and 
neutron-diffraction data. We use a similar representation as used in the 
graph for pair distribution function before: different symbols corresponds to 
the experimental data and different lines corresponds to theoretical calcu- 
lations. The theoretical structure factor curves in Figj2] agree well with the 
experimental data except in the region of diffraction maximum at 2A~^ Our 
results using the HNC approximation match the experimental results closely. 
Here again, we note that our calculation using PY equation produces a higher 
peak value while earlier classical calculations produce a lower peak compared 
to the experimental data. As we discussed in the second paragraph of this 
section, the choice of the values of the parameters in the potential function 
is responsible for this discrepancy between results from PY equation and 
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experiment. 

The structure factor is expected to approach the Feynman value [32] for 
small k, 

*^0, (27) 

where c is the velocity of sound. This limit is shown by the dash-dot line in 
Figl2j The structure factors obtained using both the PY and HNC equations 
approach a constant as k tends to zero. The reason for this discrepancy is 
that we have restricted ourselves to a pair function that remains a constant 
for large distances. 



5. Discussion 

We proposed a calculation of quantum mechanical pair distribution func- 
tion using the Percus-Yevick and Hypernetted chain integral equations. The 
interaction potential between the ^He atoms are assumed to be given by a 
hard-sphere and an attractive square-well. There is a good qualitative agree- 
ment between the theoretically obtained values of the pair distribution func- 
tion, the liquid structure factor and the experimental data. However, both of 
these computed physical quantities have higher peak values compared to the 
experimental results and the oscillations in the values are more pronounced 
for large distances. Also, the cutoff in the pair distribution function is shifted 
by a small amount. The parameters of the interaction potential adjust the 
position and sharpness of the cutoff and the peak in the pair distribution 
function. In order to get good quantitative agreement, we need to adjust 
the parameters of the interaction potential. It would be useful to incorpo- 
rate a more realistic potential between the atoms, e.g. Aziz potential[30] and 
compare the pair distribution function, obtained using our proposed method, 
with experiments. 
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Figure 1: (Color Online) Comparison of the computed quantum mechanical 
pair distribution function with earlier classical calculations and experimen- 
tal data. The solid line represents the values obtained using Percus-Yevick 
equation, while the dashed line represents the values obtained using HNC 
equation. The diamonds are the x-ray data of Gordon et al. [27]. The 
squares and the triangles show the neutron-diffraction data of Henshaw [28] 
and Svensson et al. [2^] respectively. The dash-dash-dot (McMillan |15j). 
dash-dot (Murphy [TB] using PY equation), dotted (Murphy [TB] using HNC 
equation) and dash-dot-dot (Francis [l7j) lines represent earlier classical cal- 
culations. 
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Figure 2: (Color Online) Liquid structure factor values of "^He in comparison 
with experiment and earlier classical calculations. The solid line represents 
the values obtained using PY equation, while the dashed line represents the 
values obtained using HNC equation. Our results from HNC approximation 
match the experiment closely. The diamonds are the x-ray data of Gordon 
et al. [27]. The squares and the triangles show the neutron-diffraction data 
of Henshaw [2H] and Svensson et al. [21] respectively. The dash-dash-dot 
(McMillan [15]) and dash-dot-dot (Francis [17]) lines represent earlier classi- 
cal calculations. The dash-dot line is computed using Feynman [32] theory 
with the experimental velocity of sound (267 m/sec). 
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